Assignment #3

library(tidyverse)         # for graphing and data cleaning
library(tidymodels)        # for modeling
library(themis)            # for step functions for unbalanced data
library(doParallel)        # for parallel processing
library(stacks)            # for stacking models
library(naniar)            # for examining missing values (NAs)
library(lubridate)         # for date manipulation
library(moderndive)        # for King County housing data
library(vip)               # for variable importance plots
library(patchwork)         # for combining plots nicely
theme_set(theme_minimal()) # Lisa's favorite theme

data("lending_club")
# Data dictionary (as close as I could find): https://www.kaggle.com/wordsforthewise/lending-club/discussion/170691

When you finish the assignment, remove the # from the options chunk at the top, so that messages and warnings aren’t printed. If you are getting errors in your code, add error = TRUE so that the file knits. I would recommend not removing the # until you are completely finished.

Put it on GitHub!

From now on, GitHub should be part of your routine when doing assignments. I recommend making it part of your process anytime you are working in R, but I’ll make you show it’s part of your process for assignments.

Task: When you are finished with the assignment, post a link below to the GitHub repo for the assignment. Make sure the link goes to a spot in the repo where I can easily find this assignment. For example, if you have a website with a blog and post the assignment as a blog post, link to the post’s folder in the repo. As an example, I’ve linked to my GitHub stacking material here.

LINK to Cande’s repository: https://github.com/MCTJ1998/ADS_Assignment.git

Modeling

We’ll be using the lending_club dataset from the modeldata library, which is part of tidymodels. The data dictionary they reference doesn’t seem to exist anymore, but it seems the one on this kaggle discussion is pretty close. It might also help to read a bit about Lending Club before starting in on the exercises.

The outcome we are interested in predicting is Class. And according to the dataset’s help page, its values are “either ‘good’ (meaning that the loan was fully paid back or currently on-time) or ‘bad’ (charged off, defaulted, or 21-120 days late)”.

Tasks:

1. Explore the data, concentrating on examining distributions of variables and examining missing values.

lending_club %>% 
  count(Class)

lending_club %>% 
  select(where(is.numeric)) %>% 
  pivot_longer(cols = everything(),
               names_to = "variable", 
               values_to = "value") %>% 
  ggplot(aes(x = value)) +
  geom_histogram(bins = 30) +
  facet_wrap(vars(variable), 
             scales = "free")

lending_mod <- lending_club %>% 
  mutate(Class = as.factor(Class)) %>% 
  mutate(across(where(is.character), as.factor)) %>% 
  select(-funded_amnt,
         -verification_status,
         -annual_inc) %>% 
  add_n_miss() %>% 
  filter(n_miss_all == 0) %>% 
  select(-n_miss_all)

lending_mod %>% 
  select(where(is.factor)) %>% 
  pivot_longer(cols = everything(),
               names_to = "variable", 
               values_to = "value") %>% 
  ggplot(aes(x = value)) +
  geom_bar() +
  facet_wrap(vars(variable), 
             scales = "free", 
             nrow = 2)

2. Split the data into training and test, putting 75% in the training data. Stratify by Class (add strata =Classto theinitial_split()` function).

set.seed(494) # for reproducibility

lending_split <- initial_split(lending_club, 
                               strata= 'Class',
                               prop = .75)

lending_training <- training(lending_split)
lending_test <- testing(lending_split)

3. Set up the recipe and the pre-processing steps to build a lasso model. Some steps you should take:

• Use step_upsample() from the themis library to upsample the “bad” category so that it is 50% of the “good” category. Do this by setting over_ratio = .5.
• Use step_downsample() from the themis library to downsample the “good” category so the bads and goods are even - set under_ratio = 1. Make sure to do this step AFTER step_upsample().
  
• Make all integer variables numeric (I’d highly recommend using step_mutate_at() and using the all_numeric() helper or this will be a lot of code). This step might seem really weird right now, but we’ll want to do this for the model interpretation we’ll do in a later assignment.
  
• Think about grouping factor variables with many levels.
  
• Make categorical variables dummy variables (make sure NOT to do this to the outcome variable).
  
• Normalize quantitative variables.

Once you have that, use prep(), juice(), and count() to count the number of observations in each class. They should be equal. This dataset will be used in building the model, but the data without up and down sampling will be used in evaluation.

set.seed(456)

lasso_recipe <- recipe(Class ~ .,
                       data=lending_training) %>% 

  #Pre-processing:
  step_upsample(Class, over_ratio = 0.5) %>% #This function creates a specification that will replicate rows of a data set to make the occurrence of levels in a specific factor level equal.
  step_downsample(Class, under_ratio = 1) %>% 
  step_mutate_at(all_numeric(),
                 fn= ~as.numeric(.)) %>% 
  #step_novel(all_nominal_predictors()) %>% 
  step_nzv(-all_outcomes()) %>% 
  step_dummy(all_nominal_predictors()) %>% 
  step_normalize(all_numeric_predictors())

lasso_recipe %>% 
  prep(lending_training) %>%
  # using bake(new_data = NULL) gives same result as juice()
  # bake(new_data = NULL)
  juice() 

4. Set up the lasso model and workflow. We will tune the penalty parameter.

lasso_mod <- logistic_reg(penalty =tune()) %>%
  set_engine("glmnet") %>% 
  set_mode("classification")

lasso_wf <- workflow() %>% 
  add_recipe(lasso_recipe) %>% 
  add_model(lasso_mod)

5. Set up the model tuning for the penalty parameter. Be sure to add the control_stack_grid() for the control argument so we can use these results later when we stack. Find the accuracy and area under the roc curve for the model with the best tuning parameter. Use 5-fold cv.

Tune:

set.seed(494) #for reproducible 5-fold
lending_cv <- vfold_cv(lending_training, v = 5)

# penalty grid - changed to 10 levels
penalty_grid <- grid_regular(penalty(),
                             levels = 10)

# add ctrl_grid - assures predictions and workflows are saved
ctrl_grid <- control_stack_grid()

set.seed(494) # for reproducibility

# tune the model 
lasso_tune <- 
  lasso_wf %>% 
  tune_grid(
    resamples = lending_cv,
    grid = penalty_grid,
    control = ctrl_grid)

 #Best tuning parameter ROC

lending_best_p_roc<- 
  select_best(lasso_tune,
              metric = 'roc_auc')
lending_best_p_roc

 #Best tuning parameter Accuracy

lending_best_p_ac<- 
  select_best(lasso_tune,
              metric = 'accuracy')
lending_best_p_ac

lending_final_wf <- 
  lasso_wf %>% 
  finalize_workflow (lending_best_p_roc)

lending_final_wf

## ══ Workflow ════════════════════════════════════════════════════════════════════
## Preprocessor: Recipe
## Model: logistic_reg()
## 
## ── Preprocessor ────────────────────────────────────────────────────────────────
## 6 Recipe Steps
## 
## • step_upsample()
## • step_downsample()
## • step_mutate_at()
## • step_nzv()
## • step_dummy()
## • step_normalize()
## 
## ── Model ───────────────────────────────────────────────────────────────────────
## Logistic Regression Model Specification (classification)
## 
## Main Arguments:
##   penalty = 0.0774263682681128
## 
## Computational engine: glmnet

lending_fit <- lending_final_wf %>%
fit(data = lending_training)

6. Set up the recipe and the pre-processing steps to build a random forest model. You shouldn’t have to do as many steps. The only steps you should need to do are making all integers numeric and the up and down sampling.

set.seed(456)
rf_recipe  <- recipe(Class ~ ., 
                     data=lending_training) %>% 

  step_upsample(Class, over_ratio = 0.5) %>% 
  step_downsample(Class, under_ratio = 1) %>% 
  step_mutate_at(all_numeric(),
                 fn= ~as.numeric(.))

7. Set up the random forest model and workflow. We will tune the mtry and min_n parameters and set the number of trees, trees, to 100 (otherwise the next steps take too long).

rf_model <-  
  rand_forest(mtry = tune(), 
              min_n = tune(), 
              trees = 100) %>% 
  set_mode("classification") %>% 
  set_engine("ranger")

rf_workflow <-
 workflow() %>% 
  add_recipe(rf_recipe) %>% 
  add_model(rf_model) 

8. Set up the model tuning for both the mtry and min_n parameters. Be sure to add the control_stack_grid() for the control argument so we can use these results later when we stack. Use only 3 levels in the grid. For the mtry parameter, you need to put finalize(mtry(), lending_training %>% select(-Class)) in as an argument instead of just mtry(), where lending_training is the name of your training data. This is because the mtry() grid will otherwise have unknowns in it. This part can take a while to run.

set.seed(494)

# penalty grid - changed to 10 levels
rf_penalty_grid <- grid_regular(finalize(mtry(), lending_training %>% select(-Class)), min_n(),
                             levels = 3)

# add ctrl_grid - assures predictions and workflows are saved
rf_ctrl_grid <- control_stack_grid()

set.seed(454)
rf_tune <-
  rf_workflow %>% 
  tune_grid(
    resamples = lending_cv,
    grid = rf_penalty_grid,
    control = rf_ctrl_grid)

9. Find the best tuning parameters. What are the accuracy and area under the ROC curve for the model with those tuning parameters?

#Best tuning parameter ROC

rf_best_p_roc<- 
  select_best(rf_tune,
              metric = 'roc_auc')
rf_best_p_roc

#Best tuning parameter Accuracy

rf_best_p_ac<- 
  select_best(rf_tune,
              metric = 'accuracy')
rf_best_p_ac

• In contrast with the LASSO metrics, here we can see that the best tunning parameters differ for roc_auc and for accuracy. For roc_auc we have mtry=1 while for accuracy that is mtry=11. In terms of min_n we have 40 for roc_auc and 2 for accuracy.

10. Next, we will fit a boosted tree using xgboost. We will only tune the learn_rate parameter. I have specified the model, recipe, and workflow below already (uncomment the code - you can this by highlighting it and then in the code tab at the top, choose comment/uncomment lines). You need to set up a grid of ten values for the tuning parameter and tune the model. Be sure to add the control_stack_grid() for the control argument so we can use these results later when we stack.

xgboost_spec <- 
   boost_tree(
     trees = 1000, 
     min_n = 5, 
     tree_depth = 2, 
     learn_rate = tune(), 
     loss_reduction = 10^-5, 
     sample_size = 1) %>% 
   set_mode("classification") %>% 
   set_engine("xgboost") 
 
 xgboost_recipe <- recipe(formula = Class ~ ., data = lending_training) %>% 
   step_upsample(Class, over_ratio = .5) %>% 
   step_downsample(Class, under_ratio = 1) %>% 
   step_mutate_at(all_numeric(), 
                 fn = ~as.numeric(.)) %>% 
   step_novel(all_nominal_predictors()) %>% 
   step_dummy(all_nominal_predictors(), one_hot = TRUE) %>% 
   step_zv(all_predictors()) 
 
 xgboost_workflow <- 
   workflow() %>% 
   add_recipe(xgboost_recipe) %>% 
   add_model(xgboost_spec) 

boost_grid <- grid_regular(learn_rate(),
                           levels = 10)
boost_grid

set.seed(494)
registerDoParallel() 
boost_tune <- tune_grid(
  xgboost_workflow, 
  lending_cv,
  grid = boost_grid,
  control = ctrl_grid
)

11. Find the best tuning parameters. What are the accuracy and area under the ROC curve for the model with those tuning parameters?

#Best tuning parameter ROC

boost_best_p_roc<- 
  select_best(boost_tune,
              metric = 'roc_auc')
boost_best_p_roc

#Best tuning parameter Accuracy

boost_best_p_ac<- 
  select_best(boost_tune,
              metric = 'accuracy')
boost_best_p_ac

12. Create a model stack with the candidate models from the previous parts of the exercise and use the blend_predictions() function to find the coefficients of the stacked model.

lending_stack <- 
 stacks() %>% 
  add_candidates(boost_tune) %>% 
  add_candidates(lasso_tune) %>% 
  add_candidates(rf_tune)

as_tibble(lending_stack)

lending_blend <-
  lending_stack %>% 
  blend_predictions()

type=prob #give predicted prbabilites

• Create a plot examining the performance metrics for the different penalty parameters to assure you have captured the best one. If not, adjust the penalty. (HINT: use the autoplot() function).

autoplot(lending_blend)

autoplot(lending_blend, type= 'weights')

- Which models are contributing most? The one that contributes the most seems to be the logistic regression, followed by the random forest.

13. Fit the final stacked model using fit_members().

lending_final_stack <- lending_blend %>% 
  fit_members()

lending_final_stack

## # A tibble: 5 × 3
##   member                     type         weight
##   <chr>                      <chr>         <dbl>
## 1 .pred_good_lasso_tune_1_09 logistic_reg  2.99 
## 2 .pred_good_rf_tune_1_7     rand_forest   1.47 
## 3 .pred_good_rf_tune_1_4     rand_forest   0.472
## 4 .pred_good_boost_tune_1_09 boost_tree    0.369
## 5 .pred_good_rf_tune_1_8     rand_forest   0.269

• Apply the model to the training data.

#lending_fit <- lending_final_stack %>%
#fit(data = lending_training)

preds_accuracy_train <-lending_training %>% 
  bind_cols(predict(lending_final_stack, 
            new_data= .)) %>% 
  bind_cols(predict(lending_final_stack, 
                    new_data = .,
                    type = "prob"))



#preds<-collect_predictions(preds_accuracy_train) 

#Confusion Matrix 
lending_mat<-preds_accuracy_train %>%
  conf_mat(Class, .pred_class) 

 
lending_mat 

##           Truth
## Prediction  bad good
##       bad     0    0
##       good  383 7009

summary(lending_mat)

-Compute the accuracy, construct a confusion matrix, and create a density plot with .pred_good on the x-axis (the probability of a response of “good”), filled by Class. Comment on what you see.

preds_accuracy_train%>%
  ggplot(aes(x = .pred_good, fill = Class))+
  geom_density(alpha = 0.5, color = NA)

Porcentaje del tiempo que predice un case como bueno vs. el porcentaje del timepo que predice un case como malo.

14. In the previous problem, you saw that although the accuracy was quite high, the true negative rate (aka sensitivity) was terrible. It’s common to see this when one of the classes has low representation. What we want to do now is investigate what happens in each of our models. Below I’ve provided code to investigate the lasso model (where lasso_tune is the name of my tuning step). Do similar things for the random forest and xgboost models.

lasso_tune %>% 
  collect_predictions() %>% 
  group_by(id, penalty) %>% 
  summarize(accuracy = sum((Class == .pred_class))/n(),
            true_neg_rate = sum(Class == "bad" & .pred_class == "bad")/sum(Class == "bad"),
            true_pos_rate = sum(Class == "good" & .pred_class == "good")/sum(Class == "good")) %>% 
  group_by(penalty) %>% 
  summarize(across(accuracy:true_pos_rate, mean))

rf_tune %>% 
  collect_predictions() %>% 
  group_by(id) %>% 
  summarize(accuracy = sum((Class == .pred_class))/n(),
            true_neg_rate = sum(Class == "bad" & .pred_class == "bad")/sum(Class == "bad"),
            true_pos_rate = sum(Class == "good" & .pred_class == "good")/sum(Class == "good")) %>% 
  summarize(across(accuracy:true_pos_rate, mean))

boost_tune %>% 
  collect_predictions() %>% 
  group_by(id) %>% 
  summarize(accuracy = sum((Class == .pred_class))/n(),
            true_neg_rate = sum(Class == "bad" & .pred_class == "bad")/sum(Class == "bad"),
            true_pos_rate = sum(Class == "good" & .pred_class == "good")/sum(Class == "good")) %>% 
  summarize(across(accuracy:true_pos_rate, mean))

LASSO has 72% accuracy, lambda=7.7, true negative rate= 63%

If you’d like to have a better true negative rate, which models would you choose and how would you go about doing this in a less manual way (you don’t need to write code to do it - just describe it in words). Be sure to remove the eval=FALSE when you are finished.

Shiny app

For this week, there is no code to turn in for this part. You are just going to need to think about the steps to take.

If you are new to Shiny apps or it’s been awhile since you’ve made one, visit the Shiny links on our course Resource page. I would recommend starting with my resource because it will be the most basic.

Everyone should watch the Theming Shiny talk by Carson Sievert so you can make your app look amazing.

Tasks:

In the future, you are going to create an app that allows a user to explore how the predicted probability of a loan being paid back (or maybe just the predicted class - either “good” or “bad”) changes depending on the values of the predictor variables.

For this week, I want you to answer the following questions:

1. How can you save a model you built to use it later (like in the shiny app you’ll create)?

We will use the following code save(my_model , file = ‘StackModel.rda’)

2. For shiny apps that get published (like yours will), it’s very important to have ALL the libraries that are used within the app loaded. If we were going to use the stacked model, which libraries do you think we’d need to load in our app?

library(tidyverse), library(tidymodels), library(themis), library(doParallel), library(stacks), library(naniar), library(lubridate), library(moderndive), library(vip), library(patchwork), theme_set(theme_minimal())

3. You’ll want the user to be able to choose values for each variable in the model. How will you come up with the values they can choose for quantitative and categorical data? Give one example for each, either using code or in words.

For each predictor I would set a number or text input where people can choose their values.

ui <- fluidPage( textInput( “verification_status”, “Status”, value = "“, placeholder =”Choose your status"),

Then I would add a submit button.

I would do the same thing but with numeric inputs. Maybe add a square where users can type the number they want for continuous variables.

4. You will need to populate each variable with an initial value. Which value will you choose? Is there a nice way to do this programatically (ie. with code)?

I will work on this with Lisa during office hours because I genuinely don’t know.

Function Friday problems

I will link to these separately. They will be posted by Tuesday.

Coded Bias

We will be watching some of the Coded Bias film together on Thursday. It is streaming on Netflix. Write a short reflection. If you want some prompts, reflect on: What part of the film impacted you the most? Was there a part that surprised you and why? What emotions did you experience while watching?

After forcing my housemates to watch the documentary with me over dinner, I think what impacted us the most was imagining the possibility of living in a state of full surveillance. We talked about this, and even though we know that we are being tracked and spied over the consumption we make online or the content we watch, this is done low-key and we are able to ignore it, or just not think about it. If we imagine our faces being tracked in the streets, like it happens in London or in China, it gets a bit more real and all of the sudden it is scary. It is also very scary to feel that algorithms control our lives and the way we make decisions that affect the lives and destinies of millions of people around the world. From hiring, to credit loans, to housing opportunities, these algorithms are behind many companies that control the way we live in society and they are embedded with the ideas and ideologies of the person that created it and the data it was fed with. I love how the documentary exposes the many biases behind the algorithms and how they denounce the “algorithmic harm” as a power structure where we (normal people) cannot appeal to the decisions made that harm us. I totally agree with what was said about the fact that we need to be constantly monitoring our use of algorithms. The Big Brother Watch group in London attempts to do so, but there should be a regulatory body controlling how big companies use algorithms in their activities, as well. The systemic biases that already exist in our world, get translated to the technology governing our lives and that is something really worrisome to me.

REMEMBER TO ADD YOUR GITHUB LINK AT THE TOP OF THE PAGE AND UNCOMMENT THE knitr OPTIONS.